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Abstract 

Context. Radiation from accretion discs in cataclysmic variable stars (CVs) provides fundamental information about the properties of 
these close binary systems and about the physics of accretion in general. Of particular interest are dwarf-nova outburst cycles during 
which variations of the disc properties allow a detailed study of the physical processes in accretion flows. 

Aims. The detailed diagnostics of accretion disc structure can be achieved by including in its description all the relevant heating and 
cooling physical mechanism, in particular the convective energy transport that, although dominant at temperatures < 10"* K, is usually 
not taken into account when calculating spectra of accretion discs. The disc's self-consistently calculated structure and emission allow 
testing models of dwarf-nova outbursts and accretion-disc models in general. 

Methods. We constructed a radiative transfer code coupled with a code determining the disc's hydrostatic vertical structure. 
Results. We have obtained for the first time model spectra of cold, convective accretion discs. As expected, these spectra are mostly 
flat in the optical wavelengths with no contribution from the UV, which in quiescence must be emitted by the white dwarf. The disc 
structures obtained with our radiative-transfer code compare well with the solutions of equations used to describe the dwarf-nova 
outburst cycle according to the thermal-viscous disc instability model thus allowing the two to be combined. For high-temperature 
radiative discs our spectra are compatible with models obtained with Hubeny's code TLUSTY. 

Conclusions. Our code allows calculating the spectral evolution of dwarf nova stars through their whole outburst cycle, providing a 
new tool for testing models of accretion discs in cataclysmic variables. We show that convection plays an important role in determining 
the vertical disc structure and substantially affects emitted spectra when, as often the case, it is effective at optical depths t ~ 1. The 
emergent spectrum is independent of the parameters of the convection model. We confirm that, as required by the disc instability 
model, quiescent discs in dwarf novae must be optically thick in their outer regions. In general, no emission lines are present in the 
absence of external irradiation. 

Key words. accretion,accretion discs - Stars: dwarf novae 



1. Introduction 

In most cataclysmic variable stars (CVs) matter lost by the 
Roche-lobe filling low-mass secondary star forms an accre- 
tion disc around the white-dwarf primary. In bright systems 
such as nova-like stars and dwarf novae in outburst, the disc 
is the dominant source of luminosity. In quiescent dwarf no- 
. !^ vae, although rather dim, the disc emission provides crucial 
K> , information about the physics of accretion discs, in partic- 
5_] ■ ular about the mechanisms driving the instabilities that are 
d ' responsible for the dwarf-nova outbursts. It is generally ac- 
cepted (see, however, Schreiber & Lasota| l2007h that proper- 
ties of the dwarf-nova outburst cycles are reproduced well by 
the th ermal-viscous disc instability model (DIM; see iLasotal 
1200 iL for a review), but in reality the description of some 
phases of theses cycles is far from satisfactory. This is espe- 
cially true of the quiescence state. Until now models of disc 
emission spectra have been c alculated consistently only for 
hot, bright accretion discs (see IWade & Hubenvl[T99 8l and ref- 
erences therein). In such a scheme the vertical structure of 
an accretion ring is calculated using the program TLUSDISK 
[Hubeny ( 1990, 1991), whereas the spectra are sy nthesized by the 
code S YNSPEC dHubenv. Lanz. & Jeff ervl 1991 . A different ap- 
proach to radiative transfer in accretion discs was presented by 
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[Sh aviv & Wehrsd ([1991 1. hereafter SW), but also in this case it 
was applied only to hot d iscs. An at t empt to apply the SW code 
to cold quiescent discs in lldan et~al ] (fT999h was not conclusive. 

The main obstacle in solving the vertical structure of cold 
accretion discs is the presence of convection, which was diffi- 
cult to include in radiat ive transfer codes (Hube ny, private com- 
munication); models of lWade & Hubenvl (Il998h are restricted to 
temperatures for which convectio n is unimportant. In their pio- 
neering work Pringle et alJ (Il986l) used vertically integrated disc 
structures and Kurucz's stellar spectra to test the DIM versus the 
model attributing dwarf-nova outbursts to an incr ease of mass- 
transfe r rate fr om the secondary . More recently iKromer et al.l 
(l2007l see also iKellermann et alJ (l2005h ^ have calculated non- 
LTE accretion-disc spectra throughout the outburst cycle of the 
dwarf nova SS Cyg, taking into account white-dwarf irradiation. 
They found that in order to obtain agreement with observations 
of this system in quiescence one has to assume a radial disc that 
is optically thin in its outer regions. Such a structure, however, is 
inconsistent with the DIM. The AcDc NLTE code by Nagel et alj 
62004.) used to calculate these spectra does not take into account 
the (often dominant) vertical convective energy transport. 

Of related interest are radiative-transfer codes describ- 
ing emission from cold pro to-stellar accreti on discs such as 
Hug elmever et all (l2009l) and lMin et al] ( l2009l) . Neither of these 
codes includes convection. 
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There is therefore a need for a radiative transfer code con- 
sistent with the DIM. This is the main motivation behind the 
present work. In particular we have o btained a schem e that is 
totally consistent with the code of iHameurv et al.l (il998L here- 
after HMDLH), which allows following the dwarf-nova spectra 
during various phases of the outburst with special emphasis put 
on quiescence. One should stress that the principal aim of our 
project is testing the DIM, not trying to reproduce the observed 
spectra of dwarf-nova stars at all costs. In particular our code 
should allow calculating various delays in rising to (and decay- 
ing from) outburst between light at different wavelengths (com- 
monly known a the "UV-delay") and determining what proper- 
ties of the quiescent dwarf-nova disc can be accounted for by the 
DIM). 

The HMDLH version of the DIM uses a 1 + 1 D scheme in 
which the time-dependent radial evolution equations use as in- 
put a pre-calculated grid of hydrostatic vertical structures. These 
structures are calculated using the standard equations of stellar 
structure (or the grey-atmosphere approximation in the optically 
thin case). The angular-momentum viscous transport mechanism 
is described by the ff-ansatz of Shakura & Sunvaev ( 1973). For 
a given a, M/R^ and effective temperature T^ff there exists a 
unique solution describing the disc vertical structure. Such so- 
lutions are very handy for solving the time-dependent equations 
of the disc evolution but do not produce realistic emission spec- 
tra. These can be calculated, however, by using the same input 
parameters (a, M/R^ and Teft) in a radiative-transfer code on 
the condition that the vertical structures calculated by the two 
methods are the same up to the photosphere. In this way one can 
reproduce spectra of the whole cycle of a dwarf-nova outburst. 

The outline of the article is as follows. In Section |2] we 
present the model used, with special stress put on marking the 
differences with the original SW code on which it is based. We 
discuss there also in details the opacities used. Section[3]is deal- 
ing with the solutions obtained for various types of physical 
set-ups. In this section we compare our vertical structure solu- 
tions with those obtained with the HMDLH code for a hot sta- 
tionary and cold non-stationary discs. These two models rep- 
resent the main phases of the dwarf-nova outburst cycle respec- 
tively: the outburst and quiescence. We compare the correspond- 
ing thermal-equilibria ("S-curves") in 13.11 The importance of 
convection in determining the vertical disc structure is discussed 
in detail in section 13.21 In section |4] we show our hot whole- 
disc solutions as examples of spectra emitted during the decay 
from maximum of a dwarf-nova outburst. We also show that, de- 
spite differences in treatment of vertical viscosity stratification 
our solutions cornpare v ery well with the U V spectra obtained by 
IWade & Hube"nvl (Il998 i Then in section |5] we present and dis- 
cuss the structure and spectrum of cold quiescent accretion disc. 
Finally we conclude the paper in Section|6]by shortly discussing 
future developments and necessary improvements. A report re- 
lating the progress of our work on the radiat ive-transfer cod e 
was published in the conference proceedings I dan et al.l (l2008h . 
Several solutions presented in this report were only preliminary 
and have been substantially improved, some results we found to 
be incorrect. We address these points in the present article. 

2. The physical model and basic assumptions 

Our model is based on the work o fUE aviv & Wehrsel([T99Th .The 
main improvements consist in using (modern) line opacities and 
in including the convective energy transport. As far as we know 
our code is the only working accretion-disc radiative-transfer 
programme which consistently includes convection. The last im- 



provement allows combining the radiative transfer programme 
with the DIM code of HMDLH where the time-dependent evo- 
lution of dwarf novae is calculated. 



2.1. Assumptions 

We describe an accretion disc as composed of concentric rings 
orbiting the central gravitating body of mass M at a Keplerian 
angular speed 



The disc is assumed to be geometrically thin, i.e. 

zo « R, 



(1) 



(2) 



where zq is the disc height. This assumption allows decompos- 
ing the disc equations into their vertical and radial components. 
In this paper we will solve only the vertical accretion-disc equa- 
tions and obtain the corresponding emission spectra. However, 
since the radiation of a steady-state accretion disc can be con- 
sidered as a sum of radiation from individual rings our results 
can directly applied to such configurations. The same is true 
for quasi-stationary phases of dwarf-nova outbursts. Also radi- 
ation from non-steady quiescent discs of dwarf-nova stars can 
be treated as being the sum of the emission of individual rings. 
This is true as long as radial gradients are much smaller than the 
vertical ones. During the outburst of dwarf nova, when the tem- 
perature and density fronts are propagating through the disc this 
assumption is no longer generally valid. However, since radial 
gradients are comparable to vertical gradients only in a narrow 
zone close to the fronts the error on integrated spectra cannot be 
very large 

The disc is assumed to be in a vertical hydrostatic equilib- 
rium and in LTE. 



2.2. Equations 

The vertical structure of a ring at a distance R from the centre is 
found by solving the following equations: 



- the hydrostatic equilibrium equation 
dP 



dz 



where P is the total (gas + radiation) pressure; 
- the energy balance equation 



dz \ dz /„d I dz 



(3) 



(4) 



where is the energy fl ux in the vertical (z) directi on. We 
use the a prescription of IShakura & SunvaevI ( Il973h to de- 
scribe the viscous energy generation 



(5) 



Equation (|5]i defines therefore the stratification of the viscous 
energy generation; 
- the radiative energy flux is given by 



f 

Jo 



(J(z,R, A) - B(T(z,R),A))K(A)dA = Q,,-{dF/dz)con. , (6) 



where J is the mean intensity, B the Planck function and k(A) 
is the monochromatic mass absorption coefficient; 
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- the radiative transfer equation is 



In the case of a stationary disc 



n-^I = K{S - I), 



(7) 



where / is the intensity, S the source function and n is a unit 
vector in the ray direction. This equation is solved in the 
two-stream approximation (see 12.41 ). 

As in HMDL H, the convective energy flux is calculated by 
the method of lPaczvnskil(ll969l) . Whenever the radiative gra- 
dient 



j dlnT 
[dlnP 



(8) 



is superadiabatic, the temperature gradient of the structure V 
is convective (V = Vconv)- The convective gradient is calcu- 
lated in the mixing length approximation, with the mixing 
length taken as //mi = a„^\Hp, where Hp is the pressure scale 
height: 



Hp = 



(9) 



which ensures that Hp is smaller than the vertical scale 
height of the disc. The convective gradient is found from the 
relation: 



Vad + (V„d-V,dmi' + A), 



(10) 



where Vaj is the adiabatic gradient, and Y the solution of the 
cubic equation: 



'ml 



43+T^i 



y3 + yy2 + y^Y-V = 0, 



(11) 



where Tmi = fp//mi is the optical depth of the convective 
eddies. The coefficient V is given by: 



V 



-2 



3Tml 



siKxP'C'p { d\np \ 



SUo-^T^'Hp [dlnTj 



X (Vrad - Vad). 



As in HMDLH in our model we take a. 



ml 



1.5. 



Here we should mention a caveat. The description of con- 
vection used in the present version of the code is not nec- 
essarily the correct one as it is not even clear that mixing- 
length approximation can be applied to turbulent accretion 
discs. However, in the context of the "a-disc paradigm" this 
is probably the best that can be done. It is an efficient method 
of describing an accretion disc in which most of the verti- 
cal energy transport is not radiative. Finally, as mentioned 
above, this treatment of convection is used in the HMDLH 
code which is the main reason for using it also in the present 
programme. If this approximation is applicable to accretion 
discs, there are no indications as to what the value of a^^i 
sho uld be. In solar model s values between 1 and 2 are used 
(see iHameurv et al]l 19981 and reference therein); in 13.21 we 
discuss the effect of various values ami on the vertical struc- 
ture. 

The outer boundary condition is 



Dir) 



-f 

Jo 



eff- 



(13) 



M = InR 



pVr 

Jo 



dz - constant. 



where Vr is the radial flow velocity one has 



(14) 



(15) 



but our solutions apply to any distribution of accretion rate 
M(r) in a geometrically thin accretion disc. 

2.3. Input (opacities and EOS) 

Several features of the original SW code have been improved but 
the use of modern opacities and especially of the line opacities 
is a major modification we present here in some detail. 

There are several ways of calculat ing opacities. One can ei- 
ther use co des such as P HOENIX ( Hauschildt et all 1 1 99 9*) or 
"Atlas 12" (lKuruczlll992h . or one can use the OPAL database. 
Our updated version of the SW code can use the opacities from 
'Atlas 12". However, since calculating the opacities is time- 
consuming we decided to use tabulated data from the Opacity 
Project (OP) database. To tabulate the data we used the data from 
the Opacity Project - 0PCD3.3 that was downloaded from the 
website at the Centre des Donnes de Strasbourg (CDS). For a 
given abundance mixture the subroutines mx. f, mixv.f, opfit. f 
and mixz-f allowed us (with some modifications) getting the tab- 
ulated opacities for the wavelengths chosen as our basic grid. For 
the emerging spectra from every ring we use 10,000 wavelengths 
starting from lOOA up to lO^A. The grid interval varies between 
AA = 0.25A for the most important zone such as the range be- 
tween 850 - 2300A and 1 - 2 A used for wavelengths shorter than 
800A or above 2300A. For wavelengths longer than 5000A, AA 
varies logarithmically starting at 2A and getting up to 200A for 
the very long wavelengths. All data were tabulated using a grid 
for the temperature T and electron density A^^. The indices ite 
and jne are defined by 



(12) ite - 4-0xlog(T), jne - 4-x\og(Ne), A(/fe) = A(;'ne) = 2 (16) 



The mesh in the tables is calculated for 140 < ite < 320 
therefore the lowest temperature that can be used in the code is 
3160K. Let us stress that the opacities in the present work do 
not take into account the global line broadening due effects su ch 
as velocities, expansion opacity (see lShaviv & Wehrsell2005]) or 
microscopic effects such as microturbulence. Consequently, the 
lines provide only an indication of the presence of a given ion 
and cannot be used, so far, for abundance determinations etc. 
These broadenings will be included in the code in the near future. 

The reason for allowing the code to work either with OP tab- 
ulated data or with 'Atlas 12" was to test the tabulated opacities. 
In Figure [T] we show the comparison of opacities as function of 
wavelength calculated for T = 10''K and p = 2 x 10"^ g cm"^, 
and T - 8 x lO-'K and p - 6 x 10"'"g cm""* using either "Atlas 
12" or the OP. The 'Atlas 12" continuum opacities used were 
for a solar mixture with the addition of only H and He lines. 
The opacities from OP are also for a mixture of solar abundance 
(which includes the lines and the continuum). The result of test 
shows that tabulated OP opacities are reliable. 

The two stream approximation contains the scattering coeffi- 
cient. Following the good agreement obtained between "Atlas 
12" and OP, the scattering coefficient were calculated from 
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Wavelength [A°] 




Wavelength [A ] 



Figure 1. Comparison between the opacities as function of 
wavelength obtained from "Atlas 12" of Kurucz and those ob- 
tain from OP for T = lO^K and p = 2 x 10 ^^gcm^^ (upper 



figure) and T = 8 x 10 K and p = 6 x 10" g cm""^ (lower fig- 
ure). The agreement between the two is more than satisfactory. 
(OP opacities are for solar abundance mixture while the Atlas 
ones contain just continuum plus H and He hnes.) 



"Atlas 12" taking into account scattering from H, He and elec- 
trons. One should stress, however, that scattering is not important 
in CV discs. 

The Rosseland mean opacities are taken from OPAL ta- 
bles for solar composition. While solving the radiative trans- 
fer equation we checked for the consistency between the 
Rosseland mean opacities obtained from the solution of the 
radiative transfer equation and the values from the OPAL ta- 
ble. The equation of state is inte r polate d from the tables of 
iFontaine & Graboske & Van Horn I (Il977h . The computer pro- 
gramme presented in this article is able to produce accretion disc 
models for temperatures from 3160 to 100 000 kQ. Consistently 
with the lower temperature limit we not take into account molec- 
ular processes and opacities. Rayleigh scattering is also ne- 
glected. Most of the reported quiescent dwarf-nova disc effective 
(colour) temperatures (> 3000^) justify this assumptions. 



2.4. The Solution Method 

The basic idea behind the SW original code was to couple the 
hydrostatic and the radiative transfer equations. To achieve this 
the programme iterates the two equations at the same time. The 
iteration is repeated until a structure is found to which both the 
solutions of the hydrostatic and the radiative transfer equations 
converge. 



The temperature profile is determined from Eq. (|6]l. The tem- 
perature T(z,R) in this equation is solved by expressing J in 
terms of B u sing the radiative transfer equation (for a d etailed de- 
scription see lWehrse ll98ll:lKalkofen & Wehrselll984 '). As men- 
tioned above we assumed that each radial ring is independent 
of its neighbours rings so we consider a vertical radiative flux 
only. The radiative transfer equation is solved in the two stream 
approximation. represent the specific intensity in the outward 
direction and I is the specific intensity in the inward direction. 
The radiative equation is therefore written as: 



dz 



= -{K{A) + a-{A))I^{R,z,X) 

+(t{X)J{R, z. A) + k(A)B(T(R, z). A) 



(17) 



' However, we do not attempt to produce stationary accretion discs 
spanning such a range of temperatures since such discs do not exist. 



where cr{A) is the scattering coefficient. The number of grid 
points is the same in the radiative transfer and the hydrostatic 
calculations. The boundary conditions are: (1) no incident flux 
(/Tj = 0), (2) at the symmetry plane the upward and downward 
specific intensities are equal: /J^ - Ij^^ where j designates the 
ring number 

As is now well known, when solving the vertical-structure 
equations of an accretion disc the main problem is that contrary 
to the case of stellar atmospheres the exact height of the pho- 
tosphere is unknown a-priori. Because of the gravity increasing 
with height, a change in the position of the photosphere is fol- 
lowed by a change of the gravitational acceleration which leads 
to a change in the entire vertical structure of the disc. One has to 
remember that the upper part of the disc is optically thin and (de- 
pending on the assumptions) may be the site of a non negligible 
energy production. However, since this optically thin region is 
a poor radiator and absorber, even a small energy production in 
this region may have large effect on the final temperature struc- 
ture. 

In view of the above situation we approached the problem 
in the following way: we guessed first an initial temperature- 
optical depth relation and assumed an initial height zq. The op- 
tical depth (measured vertically from outside the disc towards 
the midplane) at the disc's uppermost point was chosen to be 
between 10"^ - 10"^, and the validity of this assumption was 
checked once the final structure was obtained. Using the T(r) re- 
lation we can integrate the hydrostatic equation from zo down to 
z = 0. One should emphasize that t in the current paper is calcu- 
lated using the Rosseland mean opacity and not the optical depth 
at 5000A as in SW. When the height zo is the correct height the 
total energy generation in the ring calculated from the energy 
equation is equal to the total luminosity of that particular ring. If 
this were not the case we were altering the value of zo, keeping 
the T(r) relation fixed, until the energy condition was satisfied. 
Once this iteration converged the hydrostatic model was consis- 
tent with the energy generation but not yet with the radiation 
field. 

We then turned to solving the radiative transfer equation us- 
ing the disc height, the pressure and the density structure ob- 
tained from the hydrostatic iteration. This structure was kept 
fixed during the iteration for the radiative field. The first step 
in solving the radiative transfer equation is to calculate the con- 
vection. For each grid point we calculated the convective flux. If 
the convective flux was dominant (more than 80% of the viscous 
flux) - the temperature was determined by the adiabatic gradi- 
ent. In the opposite case the convective flux was subtracted from 
the total (viscous) flux and the solution of the radiative transfer 
equation for the remaining flux, yielded a new T(t) relation. The 
temperature profile was assigned new values for the old values of 
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Start 



Assume T(t) 

1 

Guess Zo, the height of the photosphere 

i 

Solve the hydrostatic eq. P and p 




Guess another Zo 



Calculate gradient and convection fluxes 



i 

Solve the R.T. eq. for {Fyisc — Fconv) 
and energy eq. T'new(T) 




Figure 2. The iteration scheme. 



12000 




Optical Depth 

Figure 3. Variation of the temperature-optical depth profile dur- 
ing the iteration to the converged solution (M^d = 0.6, Mq, 
R = S/^wd, a = 0.03 Teff = 5000K). The inner part of the struc- 
ture is convectively dominated up to t = 3. In this case the so- 
lution was obtained after seven iterations. In the last iteration 
the temperature profiles obtained from vertical equilibrium and 
radiation field differ by less than 0.5 % 



optical depth. If the new T(t) profile from the radiative iteration 
agreed (to a chosen degree of accuracy; usually 0.5%) with the 
old T(t) relation from the hydrostatic iteration, the hydrostatic 
structure was consistent with the radiation field. If not, we iter- 
ated again the hydrostatic equation using the new T{t) relation. 

In Figure [3] we show how the temperature profile changes 
during the iteration for a cold accretion disc. In this case the con- 
vergence was obtained after 7 iterations. One has to remember 
that since the location of the photosphere is unknown at the be- 
ginning of the iteration we must iterate for it so as to satisfy the 
energy conservation condition. The initial temperature optical- 
depth profile is a priori unknown and we guess a certain T(t) 
law. This T(t) might change dramatically after the first iteration 
but afterward quickly converges to the correct solution. 



To summarize we have four basic iteration loops in the code 
(see Fig.|2]i: 

1 . Calculating the hydrostatic structure using a given T(t) rela- 
tion. The output results are disc height, pressure and density 
structure. 

2. Calculating the convective flux using the structure of the hy- 
drostatic structure and subtracting it from the total flux. 

3. Solving the radiative transfer for a given disc height, P and 
p. The output result is a new T{t) relation. 

4. Repeating the previous two iteration until the T{t) relation 
from the hydrostatic part is in agreement with the T{t) rela- 
tion from the radiative transfer equation. 

Only after those iterations have converged the flux emerging 
from the disc is calculated using the radiative transfer equa- 
tion. The vertical structure of the ring itself is divided into 100 
z points (100 layers) and we use 10000 wavelengths to calculate 
the radiative transfer equation and the emerging spectrum. 

Once lines opacities are included in the code the calculation 
of the spectrum must be done with a much finer grid. The con- 
tribution of the lines to the total energy budget in the disc when 
solving the radiative transfer equation is negligible. Therefore 
100 grid points can be used to obtain the convergence of the hy- 
drostatic and radiative transfer equations. However, for the cal- 
culation of the emergent spectrum we interpolate between the 
100 vertical grid points to up to 500 point s in order to ge t bet- 
ter accuracy. In a recent progress report dldan et al.l 1200 8) we 
announced the presence of emission lines in the spectrum, un- 
fortunately this was an artefact of using only 100 grid points in 
the calculation. Once a finer grid was used, the alleged lines dis- 
appeared. This just indicates the importance of the accuracy in 
calculations of the spectra. 

The entire disc is divided into a series of concentric rings and 
their width is determined according to their distance from the 
WD. All rings which are further away than 2R„d are assumed 
to have a width of IR^a, below 2i?wd the width of the ring is 
taken to be 0.05/?wd, which allow us a detailed calculation of the 
emission of the rings near the boundary layer. 

The above procedure is repeated for each ring. The emerging 
fluxes from the rings are then integrated in order to get the total 
flux from the entire disc according to: 

Ftot.i = 2n F{R)RdR (18) 

Since the radiative transfer equation is solved in the two- 
stream approximation the resulting spectrum does not include 
the effect of limb darkening. The angle dependence will be added 
to our code in a future paper through a SYNSPEC-type code. 
In the present article, however, we will show spectra convolved 
with an instrumental profile in the manner of the ROTI>Q pro- 
gramme. 

2.5. The disc instability model 

Since the main objective of our work is obtaining a tool allowing 
tests of the dwarf-nova Disc Instability Model (DIM; see Laso^ 
2001, for a review) we will briefly recall its main tenets that are 
relevant to the present work. 

In the DIM the dwarf-nova outbursts are attributed to a 
thermal-viscous instability appearing when the effective temper- 
atures are contained between 5000 < Teft < 7000K (these criti- 
cal temperatures depend only very weakly on the other disc pa- 
rameters). Therefore stationary discs with effective temperatures 

" http://nova.astro.umd.edu/Synspec43/synspec-frames-rotin.html 
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in this range do not exist in the Universe. When the rate at which 
a disc is fed with matter corresponds to an unstable configura- 
tion the result is an outburst cycle during which the disc oscil- 
lates between the hot and cold (quiescent) states. Outbursts are 
triggered by heating fronts propagating up (inside-out) or down 
(outside-in) the surface-density gradients. During the rise of the 
outburst the disc is hot downstream and cold upstream of the 
front. Finally, at maximum most of the disc is brought to a hot 
quasi-stationary state with the characteristic T^ff ~ R^^l'^ pro- 
file. A cooling front (always outside-in) brings the disc back to 
the cold (quiescent) state through a sequence of quasi-stationary 
configurations with decreasing accretion rates. In quiescence 
the effective temperature is everywhere lower than characteris- 
tic temperature (Eq. l22l i and is roug hly constant with ra dius (as 
conformed by observations, see e.g. Froning et al.lll999l) . On the 
other hand the surface-density and the accretion rate strongly 
increase with radius (in quiescence the disc is filling up). The 
fact that spectra of quiescent dwarf-novae cannot be fitted with 
constant accretion-rate disc models is confirmation of this pre- 
diction of the model. All the phases of the outburst cycle have 
characteristic spectral signatures. 

Locally (for a given ring), for a given set of parameters 
(white-dwarf mass, viscosity parameter) the disc equilibria form 
an S-shaped curve on the 2 - T^s plane (see Fig.|4]l. During the 
decay from maximum the surface density and effective tempera- 
ture of the hot parts of the disc decrease along the upper branch 
of the S-curve until they reach the critical value below which 
no stationary solution exists. The value of the critical surface- 
density 
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implies that the outer disc in quiescence must be optically thick; 
especially that after leaving the upper branch its surface density 
increases before it settles on the lower (cold) branch of the S- 
curve. Therefore observations and spectral models that requires 
an optically thin quiescent outer disc in dwarf-nova stars contra- 
dict the DIM (see Sect.|5]l. 

3. Single-ring calculations 

3. 1 . The S-curves 

Since we intend to calculate spectra of disc states given by the 
time dependent code of HMDLH we checked that the vertical 
structures calculated by the two codes are sufficiently close for 
the whole procedure to make sense. We compared in detail hot 
and cold disc vertical structures. 

The principal difference between the two calculated disc 
structures is their heights. In the HMDLH case the boundary 
condition is fixed at t - 2/3 whereas in ILHS it is given by 
T^Q(~ 10-2 - 10-«). 

The best way of comparing the two structure is by plotting 
the thermal equilibria obtained by the two methods. As men- 
tioned above (Sect. 12.5b . in the (2, Teff) plane they form an S- 
shaped curve. Figure |4] shows an example of reff(S) S-curves 
obtained by the HMDLH and ILHS codes for M„d - 0.6 Mq 
at a ra dius of 5j?wd (we used the white-dwarf mass-radius rela- 
tion of lNauenber jll97 2). Each point on the S-curve represents 
a thermal-equilibrium state of the accretion disc at that radius. 
The upper branch of the S-curve corresponds to hot stable solu- 
tions, while the lower branch where the temperature is below 
- 6000 /T, represents cold, stable configurations. The middle 
branch of the S-curve contains thermally unstable equilibria. 




i:(g/cm'') 



Figure 4. E - T^ff curves at R - 5R„d with Mwd = 0.6Mq, for 
a - 0.3 (left) and 0.3 (right). In solid line the curved obtain from 
HMDLH model and in "+" the results obtained with ILHS. The 
boxes represent the zones where convection dominates. 



The agreement between the results of the two codes is very 
good on both branches of the S-curve. On top of the changes in 
the opacities due to the hydrogen ionization, the disc close to the 
turning knees of the S-curve (the unstable zone) is fully convec- 
tive (see Fig.|4|. We were unable to calculate the unstable branch 
with our radiative transfer equation code. This is of not much 
importance as the equilibrium states represented by the middle 
branch are physically meaningless. However, we were also un- 
able to obtain the convection-dominated disc vertical structures 
between T^s ~ 10"^ K and the critical hot-branch temperature 
X 7000 K since not enough radiative flux was left to use in the 
transfer code. This is a more serious drawback which somewhat 
limits the code's application. We are working on trying to fix this 
problem. 

On the hot branch, for Teff ~ 10"^ K the convective flux 
transports ~ 100% of the flux. Configurations with temperatures 
down to Teff «i 4000 K are also convection dominated. For ex- 
ample atR = 5Rv,d, for 5000 K 70 % of the disc vertical extent 
is convective and convection still transports ~ 100 % of the flux. 
For effective temperature of 4000 K 35 % of the disc has a con- 
vective flux which represent 60 % of the total flux. 

The difference between the values of Smax (the maximum 
surface density on the lower branch) calculated by the two codes 
is about 2 %. 

The excellent agreement between the low-temperature solu- 
tions calculated by the two methods shows that molecular opac- 
ities do not affect the disc structure: they are taken into account 
in HMDLH, but not in our code. 

3.2. The effect of convection 

The discussion in the preceding subsection suggest that the ef- 
fects of convection are worth a closer look. 

At temperatures corresponding to hydrogen recombination 
convection begins to play an important role in the vertical energy 
transport. Therefore at these and lower temperatures it must be 
taken into account when determining the vertical structure of ac- 
cretion discs. In Fig.|5]and Fig.|6]we present the vertical structure 
of a disc ring with effective temperature of 5000 K and a - 0.03 
at a distance of 75/?wd- Convection dominates the inner (close to 
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Figure 5. The total (viscous) flux and the convective flux (dashed 
Hne) as function of the total optical depth calculated for a cold 
disc with a = 0.03, T^^ = 5000 K, Mwd = 0.6Mq, at 7? = 75Kd 
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Figure 7. Comparison between the fluxes for models with and 
without convection for a cold disc at R=75i?wd and Teff^SOOOK 
and for a - 0.03. Convection is still efficient for r > 1 which 
aff'ects the emerging spectrum. 
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Figure 6. Comparison between ttie temperature profiles for mod- 
els with and without convection for cold disc at R=75/?wd and 
reff=5000K and for a = 0.03. 

midplane) part of the disc up to t w 10. (At Teff ~ lO^K the 
convection zone reaches t ^ 2). 

Figure |6] exhibits the influence of convection on the tem- 
perature profile. Once convection is included the vertical struc- 
ture becomes cooler and denser. Therefore the disc is becoming 
optically thicker - for the model shown here the optical depth 
changes from r w 10'* without convection, to almost t « 10^, 
when the latter is included. 

Therefore effects of convection noticeably affect the emerg- 
ing spectrum. According to Figure|6]at r = 1 the temperature of a 
convection dominated disc is lOOOK cooler than a radiative disc 
with the same effective temperature. The different temperature 



structure implies therefore different emergent spectra. Figure Q 
presents the differences between the corresponding disc spectra 
obtained with and without convection. 

The calculations presented in this paper do not include disc 
irradiation by the white dwarf and boundary layer. In principle, 
by heating up the disc surface irradiation ma y reduce the effects 
of convection (see e.g. lHameurv et al1ll999l) . However, in prac- 
tice, only the innermost regions of a qui escent dwar f-nova disc 
will be strongly modified by irradiation dHameurve t al. 19991, 
the rest of the disc being only weakly modified especially that 
the optical depth increases with radius. There ir radiation mi ght 
affect the spectrum by producing emission lines (ISmaklll991h . 

As stated in the previous section we choose = 1 .5 in or- 
der to compare our results to those of HMDLH. In order to see 
how the disc structures and the emerging spectra are affected by 
varying the mixing length parameter, we solved the equations 
assuming the two values deUmiting the range allowed by so- 
lar type stars (Demarq ue & Guentheilll99lt Guzik & LebretonI 
[l991). Figure [8] shows the vertical structures of a T^ff - 5000K 
ring at 75/?wd for a,,,/ = 1 , 1 .5 and 2. As a„i decreases convection 
becomes less efficient and the temperature therefore increases. 
Changing a,„i will affect the S-curves (discussed in the next sec- 
tion). However, since changes in the upper part of the disc due 
to changes in a,„i are minor (Figure [8] bottom), the emerging 
spectrum is hardly affected by the assumed value of the mixing- 
length parameter 

4. Full disc emission 

4.1. Solutions at high accretion rate (a = 0.3 j 

As recalled in Sect. 12.51 during the dwarf-nova outbursts ac- 
cretion discs decay from maximum through a sequence of hot 
configurations with decreasing accretion rate and outer radii 
(see e.g. Fig. 5 in Hameurv et al. 1998). Spectra of three such 
configurations are shown in Figure |9] where a - 0.3, M^d = 
0.8 Mq have been assumed. The disc is seen face-on (/ - Q") 
and the white-dwarf spectrum has not been added. The con- 
figurations correspond to three values of the accretion rate: 
M = lO'** , lO'^ , lO'^ "' g s"' (maximum temperature T^s = 
70 000 , 24 000 and 16 OOOK respectively). The outer radii corre- 
spond to the position of the cooling front, i.e. to the radius where 
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Figure 8. The vertical temperature structure for ami = 1,1-5 and 
2 (top) and as function of the optical depth (bottom). The model 
pai-ameters are: a - 0.03, M„d = 0.6, Tefj = 5000 and /?wd - 



Jeff ~ Here we chose 35, 14and3/?„d (2.45x10'", 9.8x10'' 
and 2.1 X lO** cm) . 

In aU these spectra one can identify the Balmer lines (at 
3970 A, 4 1 00 A, 4340 A, 4860 A, 6560 A) as well as the Paschen 
series and the Lyman lines. Until now we have presented un- 
convolved "raw spectra" to show the direct output of the radia- 
tive transfer calculation. Figure [10] shows the comparison be- 
tween the rest- frame spectra of an accretion disc in the UV 
part of the wavelengths and the same spectrum convolved with 
a Gaussian instrumental profile having a FWHM of 1 A sampled 
every 0.1 5 A. 

4.2. Comparison with tlie \Wade & Huben]^ accretion 
disc spectra 

IWade & Hubeiivl (l99^ hereafter WH) calculated with 
TLUSTDISK and SYNSPEC a large grid of far- and mid- 
ultraviolet spectra (850-2000 A) of the integrated light from 
steady-state accretion discs in bright cataclysmic variables. 
These data have being used extensively in the recent years for 
studying the observed accretion disc spectra. In WH model the 
disc spectra are calculated in four steps. In the first step they 
calculate the vertical structure of concentric annuli, with each 
annulus behaving as an independent plane-parallel radiating 
slab. The energy balance is enforced between radiative losses at 
the disc surface and heat generation due to viscosity. Once the 
vertical structure of the annuli is obtained, the radiative transfer 
equation is solved to compute the local, rest-frame spectrum for 
each ring of the disc. In the third stage the rest-frame intensities 
are combined to generate an integrated disc spectrum, and in 
the fourth stage the monochromatic fluxes are convolved with a 



10^ 




Wavelength (A°) 

Figure 9. The full disc spectra for three values of accretion 
rates and outer radii corresponding to the decay from outburst. 
M - 10'^ , 10'^ , 10'^-^ g s"'. The outer radius of the hot disc is 
respectively at 35 , 14 and 3 R„d- (a = 0.3, Mwd = O.8M0). 

Gaussian instrumental broadening function and then re-sampled 
uniformly in wavelength. All the published WH models are for 
disc atmospheres which are optically thick and are viewed from 
a d istance of lOOpc. 

IWade & Hubenvl (fT998l) also studied the effect of different 
instrumental resolutions (Gaussian FWHM = 0.1, 1.0 and 3.5 A) 
of the UV spectrum. The FWHM=0.1A corresponds to the res- 
olution of the FUSE spectrometers whereas 3. 5 A corresponds to 
the HUT spectrometer. According to Wade & Hubenv ( 1998) the 
differences between the 0.1 A and lA were barely noticeable, 
and they applied a FWHM of lA for their calculations. Having 
confirmed their results we also have chosen instrumental broad- 
ening with FWHM=1 A. 

Fig urefTTIshows a comparison between the I Wade & Hubenvl 
(Il998h accretion discs spectra (model cc and model y) and 
the spectra obtained with the present ILHS programme. The 
model parameters are: Model cc - M„d = 0.8Mq and M = 
3.16 X lO^'^Moyr"' and for model y - M^d = 1.21Mo and 
M - lO^^^Moyr"'. We used the same number of rings. 

The difficulty with comparing the two schemes lies in the 
formulation of the viscosity prescription. We are using the a pre- 
scription whereas WH define the viscosity coefficient through 
the value of Reynolds number of the flo w, the vertical stratifica - 
tion being controlled by a parameter ^ dKnz & Hubenv iri986l) . 
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Figure 10. Rest frame spectra (black) and spectra shown after 
convolution (red) with a Gaussian instrumental profile (see text 
for details). The parameters are the same as in Fig.|9]for/?out - 
35R„d- 

The relation between the two parameterizatio ns of disc viscos- 
ity is rather complicated a nd non-local (see H ubenv & Hubenvl 
Il998t iKromeret al. ] l2007h . The relation betw een the viscosity 
param eter a and the Reynolds number Re of iKnz & Hubenvl 
can be written as 



(20) 



where vk is the Keplerian speed and c, the speed of sound. 
The parameter f represents the vertical stratification of viscos- 
ity through: 



viz) = (vo|)' 



(21) 



to be compared with the ansatz v(z) = a P(z) used in our scheme. 

We estimated that our a - 0.3 should be roughly equivalent 
to the Re = 5000 and ^ = 2/3 of WH but the correspondence be- 
tween the two diff'erent descriptions of the viscosity mechanisms 
is not exact. Despite of that the two spectra are very similar and 
the small differences are in practice difficult to pinpoi n t. The 
comparison of the two methods presented in lldan et a n (120081) 
showed a rather serious difference around 1300A: a clear ab- 
sorption feature seen in Wade & Hubeny ( 1998) was missing in 
ILHS. After a careful study of both the code and the OP database 
we arrived to the conclusion that there is an error in the opacity 
data (silicon lines misplaced). This has been confirmed and cor- 
rected by Delahaye (private communication). Spectra in Figure 
[TT]have been calculated with the corrected opacity tables. Now 
both spectra show an absorption feature at 1300A. The remain- 
ing differences most probably are due to the differences in the 
parametrization (and stratification) of viscosity. 

5. Cold, quiescent dwarf-nova disc 

According to the disc instability model, the disc in the quiescent 
phase of the dwarf-nova outburst is characterized by a low a and 
a low (quasi-constant) effective temperature. The effective tem- 
perature in the entire disc must be lower than the critical value 
(see Lasota 2001) 




goo 1000 1100 1200 ISOO 1400 1500 I6OO 1700 ISOO 1900 2000 
Wavelength [A°] 
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Figure 11. Comparis on with spectra calculated by 
(IWade&Hubenvl ll998h. Model cc (upper panel) and model y 
(lower panel) are compared with IL HS spectra calculated fo r 
a — 0.3 and the same M„d and M as in Wa de & Hubenvl (Il998h . 
Both spectra are shown after convolution with a Gaussian 
instrumental profile. 



and the critical midplane temperature in the quiescent state must 
be lower tharQ: 



/ „ -.-0.01 

T; k 9000 Q'""- —r^, — K, 
UOiOcm/ 



(23) 



Both observations and models show that for most of the quies- 
cence duration the effective temperature profile is roughly paral- 
lel to the critical one. Therefore as examples of spectra of quies- 
cent dwarf-nova disc we calculated models with constant (with 
radius) eff'ective temperature. In such non-equilibrium discs ac- 
cretion rate is not constant but increases strongly (~ r^-^) with 
radius. 

Let us stress that low values of a (< 0.03) in quiescence are 
an essential ingredient of the DIM. This is because observations 
fix the value of the hot disc Q'(hot) to ~ 0.2 (Smak 1999) and 
Qf(hot)/Q'(cold) > 4 - 5 is required since otherwise the model 
produces no dwarf-nova outbursts. Since according to the DIM 
the values of the critical surface densities (see Eq. [T9j are rather 
high and the values of a in quiescence must be rather low, the 
model requires quiescent discs to be optically thick, especially 
in their outer regions. 

On this subject obs ervations are a t best ambiguous. For 
example, according to [littlefair et al.l (l200lh (confirmed by 



example, accoramg to ottieiair et ai.l (izuuii) {com 
ICopperwheat et al.l 12009 ^ optical observations of the quiescent 
dwarf nova IP Peg seem to suggest that the outer part of the disc 

^ Whereas the fit critical of T^ff is accurate within few percent, the 
accuracy of the fit to the critical central temperature is n o better than 
20%. This explains the difference with the fits given in iLasota et all 
i2008h . which in addition contain typing errors. The formula in Eq. Al 



(22) of this paper should read: = 8240 ^7""** M" "'* K. 
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Figure 12. The optical depth as function of the disc radii for disc 
models with a = 0.03, M^d - O.6M0 and two values of effec- 
tive temperatures reff=5000K, and reff=3500K. Above 35 /?wd 
the 3500 K disc is fully convective and its spectrum cannot be 
calculated by our method. 

in this system is optically thin. On the other hand^ Ribeiro et alJ 
(l2007h find that infrared observations imply an optically thick 
outer disc. On the other hand practically all observed dwarf- 
novae in quiescence have flat radial-temperature (< 5000 K) pro- 
files and strongly inward-decrea sing accretion rate dWood et alj 
II 98611 1989tlFroning et al.ll 19991) in excellent agreement with the 
prediction of the DIM. The presence of strong emiss ion lines is 
not an evidence for the optical thinness of the disc (Smak 1991) 
since they can be pr oduced in a chromosph ere above the disc. 
Also, as suggested bv lVrielmann et all (l2002h . the quiescent disc 
could be clumpy, with optically thin patches interspersed among 
optically thick gruels. In addition IP Pe g shows strong emission 
lines also in outburst (iMarsh & Hornd ll990). The whole sub- 
ject is therefore rather controversial ("based on current knowl- 
edge of such discs, optic ally thick regions might or might not 
be present", according to lHvnes et al.ll2005 ) so that it is clearly 
very important to be able to model correctly and consistently the 
emission spectra of quiescent dwarf novae if one wishes to solve 
this fundamental problem. 

Here we present the first spectra calculated for self- 
consistent cold disc structures with a = 0.03, Mwd = 0.6 
and a constant eff'ective temperature of 5000 and 3500K (The 
second value corre sponds to IP Peg in quiescence according to 
iRibeiro et al.ll2007h . The outer radius in the model was assumed 
to be 75iJwd- 

The 5000 K quiescent disc is everywhere optically thick 
(FigurefTSli: the optical depth in the outer regions reaching 10^. 
The disc is everywhere convection dominated. As shown above, 
the convective flux determines the vertical structure of cold ac- 
cretion discs and obviously cannot be ignored through the part 
of the disc where the optical depth is larger than 5-10. The 3500 
K disc is optically thin in its innermost regions but starting from 
~ 107?wd the optical depth r > 10 and reaches 10"^ at the outer 
edge. Figure[T3]shows the extension of the convective zones for 
cold quiescent discs. The height of the disc corresponding to op- 
tical depth of 10"^ is marked with a solid line while the height of 
the convection zone (where the convective flux is equal to half 
of the viscous flux) is marked with a dashed line. As seen in 
this figure the depth of the convection-dominated zone increases 
with radius. For the 5000 K disc it almost reaches the base of the 
photosphere for/? ~ 15R„d- The inner regions (/? < 7 R„d) of the 
3500 K disc are radiative but at 35 /?wd the convection dominated 
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Figure 13. The disc height at optical depth of 10"^ (solid line) 
and the height of the convection zone (dashed line) for two 
values of effective temperatures reff=5000K, and reff=3500K; 
(ff = 0.03, Mwd = 0.6Ma). 



zone reaches the photosphere. For larger disc radii the "leftover" 
radiative flux is too minuscule to be sensibly included into the 
radiative transfer code. In any case assuming in the framework 
of the DIM the optical thinness of the outer disc is totally ad 
hoc and does not provide any understanding of the quiescence 
according to this model. 

The total fluxes obtained for cold discs with an outer radius 
of 75 /?wd are shown in Fig. [14] and [15] The spectra in the range 
of 3000 - 10000 A are nearly flat as observed and no Balmer 
jump is present. As expected there is no significant contribution 
from the disc to the UV band. The 3500 K shows a Brackett 
jump in absorption, except for case when radiation comes from 
the innermost disc regions only; then the jump is in emission 
and weak emission lines are present in its vicinity. This reflects 
obviously the optical thinness of these regions but is unlikely to 
have observable consequences. 

As mentioned above quiescent dwarf nova spectra often 
show strong and broad hydrogen and helium emission lines. The 
spectrum supposed to represent a cold dwarf-nova disc in Figure 
[T4]has only absorption line. The temperature in the photosphere 
of the cold disc are too low to produce emissio n lines. (As men- 
tioned above, the emission lines we reported in lldan etal ] l2008i 
were an artefact of insufficient resolution of the computational 
grid). 

For low as w e were unable to produce a "corona" by the 
|Sh aviv & Wehrsd (il986) instability. Apparently heating was al- 
ways compensated by cooling even high in the disc atmosphere. 
Therefore it seems that to obtain emission lines one should take 
into acco unt the effect of irradiation by the w hite dwarf (as sug- 
gested bv lSm^[T99l . lKromer et al.1 illoOTl) have included disc 
irradiation in this context. They calculated spectra of an accre- 
tion disc in quiescence and in outburst attempting to reproduce 
observations of the dwarf-nova SS Cyg. Their model reproduces 
the hydrogen Balmer emission lines but this interesting result 
was obtained by assuming that the disc outer ring is optically 
thin. As reminded above such an assumption is in contradiction 
with the DIM results and, as well known, imphe s a valu e of the 
viscosity parameter a > I. Indeed. iKromer et alJ (l2007h assume 
in the outer ring a valu e of Re ten time lower than in the rest of 
the disc. As shown by Idan et all (11999 '). in a quiescent disc the 
requirement of optical thinness inevitably leads to high values of 
the viscosity parameter. Because of these arbitrary choices it is 
unclear what is the meaning of their results. 
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Figure 14. The integrated disc spectrum as function of the outer 
disc radius for cold discs with Teff = 5000 K (upper panel) and 
Jeff = 3500K (lower panel); a = 0.03, Mwd = O.6M0, R^a = 
8.5 X 10** cm. 
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lations of IWade & Hubenvl (Il998l) despite the different viscosity 
prescriptions respectively used. The a-viscosity prescription of 
our code is used in the Disc Instability Model and in general to 
characterize the structure of accretion discs. 

Although our code includes convective energy transport it 
fails to account for important effects that have been successfully 
taken into account by other authors. For example we assume 
LTE, whereas the AcDC code (Nagel et al. 2004) describes the 
NLTE structures. These different approaches are therefore only 
milestones on the road towards a complete description of ac- 
cretion disc emission. Unfortunately the main obstacle one will 
have to overcome is the persisting lack of a reliable description 
of the disc's "viscosity" stratification. 

Our code has been designed to work in synergy with the DIM 
code of HMDLH. This goal has been achieved for most of the 
range of the relevant parameters. We find, in agreement with the 
predictions of the DIM, that quiescent dwarf-nova discs must be 
optically thick. In absence of spontaneously formed disc corona 
emission lines observed in quiescent dwarf-novae have to be the 
result of irradiation by the white dwarf (and the boundary layer) 
of an optically thick disc as suggested by Smak ( 1991). 

In a subsequent article we will include the effects of limb 
darkening and produce a grid of disc spectra to be compared 
with observations of various types of cataclysmic variables and 
various phases of the dwarf-nova outburst cycle. In particular 
we will focus on the quiescence discs and the effect of the irra- 
diation by the white dwarf and the formation of emission lines. 
Combining ILHS with the time dependent HMDLH code will 
provide a tool for the study of dwarf nova outburst cycles. 

In particular it will be possible to compare observed qui- 
escent dwarf nova spectra with the model spectra of cold non- 
equilibrium accret ion discs and hot whi te dwarfs, instead of as- 
suming (as in e.g. lUrban & Sionll200^ , and references therein) 
that the disc in such a state is hot and in equilibrium. 

The structure of our code allows easy inclusion of various 
tabulated opacity data. In a future work we plan calculating spec- 
tra of helium discs and very cold discs in which molecular opac- 
ities are important. 
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Figure 15. The integrated disc spectrum for two low values of 
Teff = 5000 and 3500 (Mwd = O.6M0, a = 0.03). 



6. Conclusions 

We have constructed a numerical code solving radiative transfer 
in accretion discs for a range of effective temperatures reaching 
values as low as 3000 K. As far as we know, our code is the only 
one taking into account convective energy transport, including 
the cases when convection is dominating the energy transport. At 
high temperatures, in the case of fully radiative accretion discs, 
our spectra are in a very satisfactory agreement with the calcu- 
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